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This  document  is  the  final  report  of  a one  year  program  of  work 
concerned  with  finding  ways  to  upgrade  automated  processing  of  Defense 
Meteorological  Satellite  Program  (DMSP)  imagery  data  at  Air  Force 
Global  Weather  Central  (AFGWC).  In  the  first  of  two  states  of  work, 
surveys  were  conducted  to  identify  needs  for  imagery  processing  at 
AFGWC  and  to  find  the  best  technique  for  AFGWC  to  adopt  next  in 
upgrading  its  capability  for  automated  imagery  processing.  The 
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procedure  and  results  of  the  first  stage,  including  the  recommendation 
that  spectral  analysis  be  adopted  for  improved  feature  abstraction  in 
the  current  3DNEPH  program,  are  separately  reported.  The  procedure 
and  results  of  Stage  2,  which  focused  on  demonstrations  and  checks  of 
spectral  analysis  on  selected  DMSP  imagery  samples,  are  presented  here 
along  with  a summary  of  conclusions  and  recommendations  from  the 
contract  effort  as  a whole. 
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1.  INTRODUCTION 


Tnis  document  is  the  final  report  of  work  conducted  under 
support  of  Contract  F 19628-76-C-O 1 24  . The  contract  was  concerned 
with  finding  techniques  to  improve  automated,  realtime  processing 
of  Defense  Meteorological  Satellite  Program  (DMSP)  imagery  data 
at  Air  Force  Global  Weather  Central  ( AFGWC ) , Offutt  AFB, 
Nebraska.  The  work  was  conducted  in  two  stages.  In  the  first 
stage,  surveys  were  carried  out:  1)  to  determine  AFGWC  needs 
for  DMSP  imagery  processing;  2)  to  find  current  state  of  the  art 
techniques  in  automated  imagery  processing  relevant  to  meeting 
those  needs;  and  3)  to  recommend  one  or  more  techniques  that 
would  be  most  promising  for  AFGWC  to  adopt  next  in  its  continuing 
effort  to  upgrade  its  capabilities  in  this  critical  area  of 
meteorological  data  processing.  In  the  second  stage  the 
recommended  technique  was  demonstrated  on  selected  samples  of 
DMSP  imagery. 

The  procedures  and  results  of  Stage  1 have  been  separately 
reported  (Pickett  and  Blackman,  1976).  The  rationale,  procedures 
and  results  of  Stage  2 along  with  summary  and  conclusions  for  the 
contract  effort  as  a whole  are  the  primary  subjects  of  this 
report . 

2.  BACKGROUND 

For  background  to  the  presently  reported  work,  we  present 
first  a brief  sketch  of  the  conclusions  and  recommendations  from 
the  work  in  Stage  1 . Issues  that  shaped  our  approach  to  Stage  2 
are  discussed,  and  the  resulting  structure  and  emphasis  of  the 
work  are  described. 


; 
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2.1  Results  of  Stage  1 


The  work  in  Stage  1 resulted  first  in  two  general 
conclusions  and  related  recommendations.  The  first  conclusion  to 
be  drawn  was  that  in  those  areas,  such  as  the  forecasting  area, 
where  imagery  processing  is  entirely  visual,  needs  for  automation 
cannot  be  adequately  determined  without  conducting  systems 
analyses  and  information  flow  studies  falling  well  beyond  the 
scope  of  the  present  contract.  It  was  recommended  that  such 
studies  be  undertaken. 

The  second  conclusion  was  that  AFGWC  lacked  a fully 
adequate  long  term  strategy  for  developing  an  automated  system, 
particularly  in  terms  of  expanding  its  capabilities  to  encompass 
processing  now  done  visually.  It  was  recommended  that  such  a 
long  term  strategy  be  developed,  one  that  allows  for  the 
integration  of  visual  and  automated  processing,  probably  through 
interactive  computer  techniques,  to  provide  the  basis  for  an 
evolutionary  approach  to  full  automation. 

As  a consequence  of  these  two  conclusions,  it  was  decided 
early  in  the  contract  schedule  that  the  most  productive  effort 
would  be  one  that  focused  on  AFGWC' s present  approach  and  thrust 
in  automated  imagery  processing  as  represented  in  the  currently 
operating  3DNEPH  program.  Analyses  showed  that  there  is  a clear 
and  pressing  need  there  to  improve  imagery  analysis  at  the 
feature  abstraction  stage,  to  upgrade  capability  at  the  front  end 
of  the  system  where  those  features  that  get  passed  on  to  the 
pattern  recognition  component  get  culled  from  the  full  data 
stream.  For  several  statistical  and  practical  reasons,  spectral 
analysis  was  recommended  as  the  most  promising  technique  for 
improved  feature  abstraction. 


2.2  Issues  Snaping  the  Work  in  Stage  2 

Tne  point  of  the  second  stage  of  work  was  to  conduct  a 
feasibility  analysis  of  the  recommended  technique:  first  to 
demonstrate  its  consistency  with  AFGWC  needs,  in  other  words,  to 
show  its  usefulness  to  AFGWC;  and  second  to  show  feasibility  of 
its  Deing  interfaced  with  AFGWC  data  analysis  systems  present  and 
proposed.  How  that  goal  was  to  be  achieved  within  the  limited 
contract  time  and  resources  reserved  for  it  was  to  depend  very 
largely  on  the  particular  technique  that  was  recommended.  Our 
approach  with  respect  to  the  particular  recommendation  of  using 
spectral  analysis  for  improved  feature  abstraction  was  affected 
by  two  main  considerations. 

The  first  consideration  was  that  little  more  could  be  done 
from  an  analytical  standpoint  to  demonstrate  the  usefulness  of 
spectral  features  for  cloud  classification  than  had  already  been 
done  in  Stage  1.  As  our  surveys  revealed,  there  are  no  extant 
mathematical/statistical  models,  or  even  an  adequate  body  of 
structural  descriptions  for  telling  what  features  of  satellite 
images  are  critical  for  distinguishing  among  cloud  classes.  Part 
of  AFGWC' s long  term  strategy  for  automation  certainly  should  be 
to  look  recurrently  for  such  aids  and  provide  for  ways  to 
integrate  them  into  the  system  as  they  become  available.  In  the 
meantime,  however,  there  is  really  no  adequate  way  to  tell  how 
useful  a particular  algorithm  for  feature  abstraction  will  be 
witnout  conducting  empirical  studies,  and  those  studies  have  to 
be  substantial  enough  to  permit  valid  generalization  of  the 
findings  to  field  operations.  Such  studies  have  to  be  made  on 
samples  of  imagery  for  every  cloud  class  that  occurs  with 
appreciable  frequency,  in  order  to  test  for  all  potentially 
significant  errors  of  classification.  And  for  every  one  of  those 
classes,  the  sample  has  to  be  sufficiently  large  and  well  drawn 


to  be  validly  representative  of  that  class.  It  is  not  just  the 
numbers  of  imagery  samples  required  that  define  the  magnitude  of 
the  effort  involved,  it  is  also  the  need  for  adequate  care  in 
determining  the  "true"  classification  of  each  image.  Various 
options  are  open  for  reducing  the  magnitude  of  such  studies,  e.g. 
by  explicitly  restricting  the  geographical  or  seasonal  domain  of 
sampling,  or  by  restricting  the  domain  of  explicitly  tested 
classifications,  but  even  a highly  pared  version  of  such  a study 
is  clearly  beyond  the  scope  of  this  contract  effort.  It  was 
decided,  therefore,  to  focus  the  work  in  Stage  2 on  the  second 
and  more  limited  aspect  of  feasibility — determining  whether 
spectral  analysis  was  compatible  with  AFGWC  data  analysis 
systems.  It  certainly  was  possible  within  the  remaining  contract 
effort  to  demonstrate  spectral  analysis  of  DMSP  imagery  taken 
from  the  AFGWC  data  analysis  system.  A demonstration  of 
compatibility  at  that  level  was  set  as  the  minimum  goal  for  Stage 
2. 

The  second  consideration  was  that  as  a result  of  our  first 
report  and  related  briefings,  AFGWC  decided  to  undertake  a 
large-scale  empirical  study  of  the  type  needed  for  an  adequate 
test  of  spectral  analysis.  It  was  decided  that  AFGL  would  assist 
in  this  endeavor  by  providing  a computer  programmer  to  do  the 
encoding  necessary  to  get  spectral  analysis  up  and  running  on  the 
AFGwC  data  analysis  system.  The  test  was  to  be  conducted  as  soon 
as  the  necessary  computer  coding,  experimental  designing  and  data 
sampling  could  be  achieved.  The  test  was  to  be  designed  to 
assess  the  impact  of  spectral  analysis  on  accuracy  of  cloud 
classification  in  the  3DNEPH  context. 

In  view  of  that  decision,  it  was  determined  that  BBN  should 
coordinate  its  activity  in  Stage  2 with  that  of  the  AFGL 
programmer.  The  intent  was  to  insure,  insofar  as  possible,  that 
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ill  computer  coding  for  bBN ' s demonstrations  and  tests  would  be 
directly  transferable  to  the  AFGWC  data  analysis  system.  Thus,  a 
smooth  connection  was  sought  between  BBN ' s work  in  Stage  2 and 
that  work  on  the  large-scale  empirical  study. 

Subsequent  to  setting  these  goals,  a decision  was  made  to 
add  an  item  of  work  to  the  contract  concerned  with  estimating 
spatial  power  spectral  characteristics  of  a range  of  pure  cloud 
types.  The  point  of  this  activity  was  that  it  might  provide 
useful  leads  in  the  search  for  discriminating  spectral  features. 
Tne  plan  was  to  arrive  at  purely  qualitative  characterizations  of 
the  power  spectra  based  on  visual  inspection  of  typical  examples 
of  pure  cloud  types  as  given  by  Conover  (1962).  Thus,  this 
activity  was  added  to  the  work  in  Stage  2. 

2.3  Structure  and  Emphasis  of  Work  in  Stage  2 

The  main  emphasis  of  work  in  Stage  2 was  on  guiding 
implementation  of  spectral  analysis  and  in  conducting  selected 
checks  and  demonstrations  on  DMSP  imagery  data. 

As  regards  the  added  item  of  work,  an  effort  was  made  to 
determine  power  spectra  for  a range  of  cloud  types,  but  it  was 
judged  that  that  was  not  feasible.  First,  the  intent  of  that 
effort  was  that  the  estimates  be  obtained  from  purely  visual 
inspection  of  available  cloud  imagery  samples  (as  given,  e.g.,  by 
Conover,  1962)  in  comparison  to  standard  test  patterns  for  which 
spectra  had  already  been  computed.  Detailed  consideration  of 
that  approach  led  to  the  conclusion  that  there  was  an 
insufficient  body  of  visually  comparable  standard  test  patterns 
in  the  literature,  for  which  spectra  had  been  computed,  to  permit 
such  purely  visual  estimates.  Second,  consideration  was  given  to 
the  possibility  of  computing  spectra  on  the  imagery  samples  given 


by  Conover  but  this  also  was  judged  to  be  not  justified,  even  by 
optical  means,  for  two  main  reasons.  First,  there  was  no  way  to 
guarantee  that  spectral  signatures  obtained  for  those  imagery 
samples  would  generalize  to  the  DMSP  images  taken  at  different 
altitudes,  by  different  sensors,  and  subjected  to  corrections  not 
common  to  the  Conover  images.  Second,  it  was  judged  that  any 
leads  on  critically  discriminating  spectral  features,  obtained 
from  comparison  of  spectra  based  on  just  the  few  samples  provided 
by  Conover,  would  be  very  suspect.  Estimates  of  error  of  the 
spectral  signature  would  have  to  be  available  to  obtain  leads 
with  any  acceptable  degree  of  reliability.  Better,  it  was 
judged,  to  wait  for  data  obtained  in  the  large  scale  empirical 
study.  Thus,  the  effort  to  estimate  spectra  was  carried  only  to 
the  point  where  it  was  judged  to  be  infeasible,  and  greater  time 
was  spent  on  the  main  line  of  activity. 


One  aspect  of  the  main  Stage  2 effort  which  required  more 
attention  than  was  originally  supposed,  was  communication  with 
AFCWC , particularly  in  planning  the  large  scale  empirical  study. 
Joint  planning  was  necessary  for  BBN  to  provide  adequate  guidance 

to  coordinate  its 
analysis,  as  far  as 


to  the  AFGL 
demonstrations 


programmer  and  for  BBN 
and  checks  of  spectral 


possible,  with  AFGWC's  preparations. 


3.  DEMONSTRATIONS  OF  SPECTRAL  ANALYSIS  ON  SELECTED  SAMPLES 
OF  DMSP  IMAGERY  DATA 


Implementation  of  spectral  analysis  in  code  compatible  with 
AFGC  data  analysis  systems  was  primarily  the  responsibility  of 
AFGL . BBN  has  contributed  advice  and  direction  to  that  work,  and 
modified  certain  of  the  reported  approaches  (Booth,  1973,  Sikula, 
1 9 7 ^ ) , particularly  as  regards  the  computation  of  wave  number 
spectra.  Setting  up,  overseeing  and  interpreting  the 
demonstration  of  spectral  analysis  on  selected  samples  of  DMSP 
imagery  data  was  BBN's  other  primary  responsibility. 

In  this  section,  we  discuss  first  the  rationale  behind 
and  the  logic  of  implementing  spectral  analysis  and  spectral 
feature-based  classification.  We  then  present  the  procedure  and 
results  of  demonstrations  and  tests  on  selected  samples  of  DMSP 
imagery  data. 


i . 1 Rationale 


Clearly,  the  immediate  point  of  implementing  and  checking 
spectral  analysis  in  this  effort  was  to  enable  AFGWC  to  proceed 
with  direct  tests  of  its  utility — perhaps  as  the  key  to  a stand 
alone  satellite  data  processor--via  the  planned  large  scale 
empirical  study.  But  the  primary  importance  of  developing  power 
spectral  analysis  within  the  3DNEPH  framework  is  not  simply  to 
Duila , immediately,  a better  processor.  More  important,  it  is  to 
assess  the  potential  value  of  going  in  this  direction  of  signal 
processing,  to  gauge  the  long  run  benefits  and  costs, 
improvements  and  limitations  of  this  next  level  of  processing 
capability.  The  operational  success  of  using  spectral  features 
for  classification  depends  foremost  on  the  presence  of 
information  in  the  recorded  data  stream;  no  transformation  of 
that  data  can  create  the  required  information  if  it  is  absent. 
For  certain  classification  tasks,  HR  data  may  provide 
sufficiently  fine  detail,  while  for  other  tasks  even  VHR  data  may 
not  have  adequate  resolution.  Similarly,  information  required 
for  classification  may  be  contained  separately  or  jointly  in  the 
visible  and  infrared  (IR)  spectral  regions,  or  it  may 
reside  in  a region  (e.g.,  microwave)  not  presently  sensed. 
Identifying  and  understanding  processor  behavior  in  cases  such  as 
these  as  they  occur  in  the  operational  scenario  will  provide 
guidance  towards  the  continuing  improvement  of  the  satellite  data 
processor . 
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$.<t  Implementation 


The  conversion  of  raw  data  to  labeled  cloud  type  will 
involve  the  computation  of  power  spectra,  conversion  of  these 
data  to  wave  number  spectra  (as  discussed  below)  and  computation 
of  a statistic  sufficient  for  classification. 


3-2.1  Two-Dimensional  Spectra  and  the  Fast-Fourier 
TransformC  FFT ) 

Tne  Discrete  Fourier  Transform  ( DFT)  in  two  dimensions  is 
defined  as 


(VV 


M 

r a [ j (Ax) , 
k=l 


-i2H  [jf  (Ax)  + kf  (Ay)] 
k (Ay) ] e x y 


(1) 


with  at~*-a  representing  Fourier  transform  pairs. 

If  we  let  N=M,  Ax  = Ay  = T/N  (with  T the  analysis  window 
size),  and  compute  fx,  fy  only  at  integer  multiples  of  T 1 , then 


a (n,m) 


N 

I 

k=l 


(2) 


which  is  the  defining  form  of  the  FFT . The  inverse  FFT  need  not 
concern  us  here  because  it  will  not  be  used  below.  Note  that  for 
n=m=0,  we  have 


i (o.o)  - Z Z a (i|  . S* 
j = l k=l  \ N N 

(3) 


p 

which  is  the  sum  over  the  data  and  therefore  N times  the  average 
value  over  the  TxT  data  window.  When  computing  wave  number 

spectra,  we  may  prefer  to  normalize  Equation  2 by  dividing  by  N2 
so  that  the  value  of  the  spectrum  at  (0,0)  will  be  the  data 
average.  Although  we  have  taken  this  option  in  our  initial 
work,  there  is  a basic  normalization  dichotomy  for  deterministic 
vs.  random  processes.  Involved  here  is  a conflict  between 
coherent  and  incoherent  summation.  Further  discussion  is 
deferred  to  the  section  treating  processing  results. 

Computing  Equation  2 using  the  FFT  algorithm  produces  an 
array  of  data  representing  the  first  quadrant  in  the  two 
dimensional  spectral  plane,  because  of  the  spectrally  periodic 
nature  of  the  FFT,  the  entire  spectral  plane  can  be  generated 
from  first  quadrant  data,  allowing  direct  computation  of  wave 
number  spectra.  Figure  1 shows  schematically  a 1 6x 1 6 array 
comprising  the  original  8x8  (upper  righthand  corner)  array  of 
data,  periodically  extended.  The  numbers  in  the  arrays  may  be 
considered  as  the  indices  (relative  frequencies)  in  a FORTRAN 
array  dimensioned  (0/7,  0/7),  with  the  comma  separating  the 

indices  implied.  That  is,  the  number  70  would  refer  to  location 
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Schematic  representation  of  FFT  output  and  its  periodic 
extension.  The  block  in  the  upper  righthand  corner 
comprising  64  elements  represents  the  64  points 
generated  by  an  FFT  of  an  8x8  array.  Each  number 
in  this  subarray  may  be  converted  to  the  indices  for 
a given  spectral  component  by  inserting  a comma 
between  the  two  digits,  e.g.,  00  is  (0,0),  43  is 
(4,3) , etc. 
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7,0,  while -the  number  57  refers  to  5,7,  etc.  Equivalently,  the 
number  70  gives  the  octal  location  (decimal  56)  of  the  element  in 
question  for  a unidimensional  storage  (0/63),  e.g.,  a(0,0)  is 
found  in  location  0,  a(2,6)  is  in  decimal  location  22,  and  a(7,7) 
is  in  decimal  location  63. 

3.2.2  Wave  Number  Spectra 

The  need  for  computing  wave  number  spectra  or  some  variant 
thereof  arises  from  the  fact  that  we  wish  to  insensitize  the 
classifier  to  the  orientation  of  a given  cloud  pattern.  In 
Sikula  (1974),  and  also  apparently  from  his  computer  program,  the 
wave  number  spectral  components  are  computed  as  the  sums  of 
spectral  amplitudes  (|a(n,m)|)  for  those  elements  falling  within 
corresponding  annular  bands  in  the  frequency  plane.  Each  band  is 
one  unit  in  width  and  symmetrically  disposed  radially  about 
integer  multiples  of  the  radial  "fundamental"  frequency  1/T 
That  is,  at  least  for  n,m  < N/2,  compute 


V + . 5 


(4) 


where  extracts  the  greatest  integer  not  exceeding  the 
expression  witnin,  and  add  |a(n,m)|  to  the  pth  wave  number. 
For  n 2 N/2,  or  m 2 N/2,  or  n,m  2 N/2,  there  is  some  question  as 
to  what  actual  procedure  was  used.  It  is  important  to  note, 
however,  that  wave  numbers  computed  from  FFT  indices  directly  are 
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not  referenced  to  the  principal  part  of  the  (aliased)  spectrum 
when  either  index  exceeds  the  Nyquist  frequency.  Consequently, 
direct  use  of  Equation  4 and  elimination  of  those  spectral 
components  for  which  p > N/2  (or  worse,  p>  N/2),  will  eliminate 
much  useful  data,  particularly  those  components  in  the  second  and 
fourth  quadrants,  as  can  be  seen  in  Figure  1 (e.g.,  a(2,5), 
a(6,3). 

To  properly  compute  wave  number  spectra  as  sums  over  annular 
regions,  we  subtract  any  given  index  from  N if  it  exceeds  N/2: 
n'  = min(n,N-l)  (5) 

m*  = min(m,N-m)  (6) 

where  min(  ) extracts  the  smaller  of  the  two  volumes  within  the 
parenthesis . 

Equation  (4)  may  then  be  applied  to  n»  and  m'  to  determine  the 
proper  total  wave  number.  Here,  we  are  merely  exploiting  the 
periodic  nature  of  the  FFT  and  representing  n.  and  m by  their 
values  within  the  principle  part  of  the  spectrum  ( lfxl<  N/2,|fy| 
£ N/2  ).  Although  we  require  this  shifted  data  to  fall  within 

the  Nyquist  region  or  principle  spectral  part,  this  requirement 
does  not  simply  limit  the  useful  range  of  p to  p<N/2,  since  it  is 
clear  that  fx  = fy  = N/2  yields  p = N . In  essence,  aliasing 
occurs  above  frequencies  given  by  the  square  region  constraint 
(|f  | » |f„|  1 N/2)  and  not  by  any  equivalent  one-dimensional 
constraint  on  p.  Thus,  we  ultimately  use  the  NxN  spectral 
components  of  the  first  quadrant,  i.e.,  all  the  FFT  data,  but 
combine  their  magnitudes  based  on  first  applying  Equations  (5) 
and  (6),  and  then  Equation  (4). 

Using  Sikula's  method  described  above,  magnitudes  of 
spectral  components  summed  within  annular  bands  become  wave 
number  spectra.  Booth  (1973)  used,  various  sums  of  normalized 
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power  spectral  components  as  his  wave  number  spectra.  The  actual 
normalization  used  by  Booth  is  not  explicitly  given.  Our 
approach,  while  similar  in  spirit  to  these  two,  differs 
numerically  because:  1)  we  compute  the  root  mean  square  (rms) 
value  of  the  amplitude  of  those  spectral  components  in  each  band 
(equivalently,  the  square  root  of  the  average  power),  and  2)  all 
components  of  the  two  dimensional  spectrum  are  utilized.  One 
problem  with  using  sums  rather  than  averages  for  wave  number 
spectra,  which  will  arise  (as  it  does  in  Sikula's  data)  when 
power  spectral  covariance  data  are  not  used  to  scale  the  spectra 
(see  below),  is  that  meaningless  fluctuations  are  incorporated  in 
the  wave  number  spectra  due  solely  to  the  manner  in  which  points 
on  a Cartesian  grid  fall  within  the  annular  spectral  bands.  A 
second  problem,  again  mitigated  by  the  use  of  covariance  data,  is 
that  summing  over  annular  bands  emphasizes  high  spatial 
frequencies  and  distorts  the  spectrum  in  a manner  akin  to 
differentiating  the  input  data.  While  there  are  arguments  on 
both  sides  of  this  issue,  we  suggest  that  averaging  rather  than 
summing  produces  spectra  of  greater  intuitive  and  diagnostic 
value.  Figures  2 and  3 serve  to  illustrate  this  point. 

j 

Figure  2 shows  wave  number  spectra  for  these  cases,  with 
data  taken  directly  from  Sikula's  report.  Case  I is  reported  as 
broken  low  level  cloudiness,  Case  II  as  small  scale  cumulus  and 
Case  III  is  labeled  wispy  cirrus.  From  a textural  standpoint, 
images  in  the  report  show  Case  II  to  contain  many  small  cloud 
elements,  Case  I to  contain  many  small  and  several  large 
elements,  and  Case  III  to  be  of  a slightly  mottled  but  mostly 
uniform  nature.  Replotting  the  data  after  dividing  each  summed 
wave  number  spectral  component  by  the  number  of  terms  that  were 
used  in  the  sum,  we  have  converted  Sikula's  data  to  average  wave 
number  spectra,  as  shown  in  Figure  3.  These  are  not  rms  values 
because  we  do  not  know  the  individual  data  values  within  each 
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Figure  2.  Three  wave  number  spectra,  after  Sikula  (1974). 


-15- 


AVERAGE  WAVE  AMPLITUDE 


I I 


Figure  3.  Three  wave  number  spectra  modified  from  Sikula’s 
data.  Total  wave  amplitude  has  been  normalized 
by  the  number  of  terms  in  each  wave  number  sum. 
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band,  but  only  their  sum.  In  addition  to  a great  increase  in 
smoothness,  all  the  curves  of  Figure  3 evidence  a general  low 
pass  characteristic  we  most  often  would  expect  to  see  in  the 
spectra  of  clouds  and  other  similar  textures.  Note  that  the 
spectra  shown  in  Figure  3 are  easily  understood  in  terms  of  the 
imagery.  That  is,  the  spectrum  of  Case  I is  consistent  with  a 
mixture  of  small  and  large  elements,  Case  II  is  consistent  with 
small  elements,  and  Case  III  indicates  little  textural  variation 
beyond  wave  number  four.  In  fact,  for  Case  III  it  appears 
plausible  that  the  spectral  constancy  above  wave  number  ten  is 
not  a direct  function  of  the  cloud  condition  but  is  indicative  of 
spectrally  white  system  noise.  Comparing  wave  number  spectra  for 
Cases  I and  II,  we  see  a strong  similarity  for  spectral 
components  above  wave  number  six  which  is  reinforced  by  a 
corresponding  similarity  in  the  small  cloud  elements  of  the  two 
images. 
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3-2.3  Classification  of  Cloud  Types 


Converting  radiometric  data  to  wave  number  spectra  is  a data 
reduction  process  which  can  also  be  thought  of  as  producing 
spectral  features  for  input  to  a classifier.  In  Pickett  and 
Blackman  (1976)we  briefly  treat  linear  discriminant  classifiers 
as  particularizations  of  the  more  general  quadratic  discriminant 
classifier  and  related  to  the  concept  of  minimum  distance. 
Considering  wave  number  spectra  as  features  and  treating  the 
class  conditional  probability  densities  for  the  feature  vectors 
as  multivariate  normal,  we  have  (see  Duda  and  Hart,  1973,  for 
example) 


- ^ (a  - y_i)t  E?1  (a  - ui)  - \ log  | E.J  + log 


(7) 


Here,  §.  is  toe  feature  vector  (wave  number  spectra  alone  or  in 
concert  with  first  order  statistical  measures),  denotes 

matrix  transpose,  gA(&)  the  discriminant  for  the  ith class,  the 
mean  vector  for  that  class,  P(c4)  the  a priori  probability  of  ith 
class  occurrence,  and  £ the  covariance  matrix  with  elements 
given  by 


*jk  - E C<*j  - V(ak  * V1 


(8) 
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The  operator  E[  ] is  the  expectation  with  respect  to  the 
ensemble.  Expanding  Equation  (7)  yields 


g±  (?) 


\ E"1  a + (I"1  ui)ta  + 6. 


(9) 


with 


‘ \ Hi 


- \ log 


| E. | + log 


P(ci) 


(10) 


The  general  form  of  Equation  (9)  indicates  quadratic  discriminant 
functions . 


In  our  prior  report,  we  note  that  if  2^  is  equal  to  a2  I (I 
the  identity  matrix),  or  even  if  2^  is  the  same  for  all  classes 
(=E,  say),  then  the  discriminant  surfaces  are  hyperplanes  and  the 
discriminants  are  linear.  To  construct  linear  discriminants, 
Booth  used  the  average  of  the  covariance  matrices  for  all  classes 
as  the  covariance  matrix  for  each  class.  However,  we  are  not  so 
much  interested  in  linearizing  the  classification  procedure  as  we 
are  in  simplifying  the  computations.  Since  we  need  to  compute  0. 

— 1 t ^ 

once  for  each  class  and  can  compute  once  and  store  (E^  p^) 

as  a new  vector  ^ for  each  class,  the  potential  computational 
complexity  resides  in  the  first  term  on  the  righthand  side  of 
Equation  9.  While  in  general,  an  NxN  matrix  will  require  N2+N 
multiplications  and  a similar  number  of  additions  (for  each 
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class),  if  is  a diagonal  matrix  we  can  scale  each  element  of  a 
by  its  class  conditional  standard  deviation  and  compute  the  norm 
of  the  scaled  vector.  This  latter  computation  will  require  only 
2N  multiplications  and  N additions  for  each  class,  which  should 
hardly  be  noticed  when  added  to  the  FFT  and  wave  number 
computations.  Assuming  or  even  forcing  a diagonal  form  for  is 
quite  reasonable,  given  the  asymptotically  independent  nature  of 
power  spectral  estimates  and  the  fact  that  the  data  window  has  no 
spectral  leakage  at  the  FFT  frequencies,  i.e.,  the  data  window 
does  not  introduce  dependence  between  different  spectral 
components.  Thus  while  it  would  be  nice  if  all  class  conditional 
covariance  matrices  are  approximately  equivalent  to  each  other, 
this  condition  is  not  required  for  rapid  classification. 
Accepting  the  covariance  matrices  as  yiey  are  rather  than 
averaging  them  (or  ignoring  them)  should  enhance  the  accuracy  of 
the  classification  procedure. 
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3-3  Processing  Cloud  Image  Data 

During  the  later  stages  of  this  program,  as  the  capability 
to  compute  wave  number  spectra  was  developed  at  AFGL,  we 
attempted  to  perform  a pilot  study  on  selected  weather  imagery. 
Unfortunately,  the  images  available  were  from,  display  tapes  and 
represented  only  HR  data.  Also,  while  three  areas  were  selected 
containing  mixed  stratus  and  strato  cumulus  in  the  first,  stratus 
or  fog,  in  the  second,  and  cumulus  and  cumulosnimbus  in  the 
third,  difficulty  in  reading  the  Cu-Cb  area  on  the  tape  left  us 
with  two  areas  containing  somewhat  similar  basic  cloud  types.  It 
was  decided,  therefore,  to  use  the  two  data  sets  as  a means  to 
test  and  debug  programs  and,  hopefully,  to  demonstrate  several  of 
the  elements  of  the  spectral  analysis-classification  scheme. 
Actual  classification  was  not  performed  because:  (1)  the  computer 
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programs,  although  straightforward,  have  not  been  written, 
( 2)climatology  data,  required  for  Bayesian  classification,  were 
not  yet  available,  and  (3)  using  data  from  mixed  and  overlapping 
cloud  types  as  a training  set  is  not  simple  nor  recommended. 

3.3.1  Average  Wave  Number  Spectra 

Figure  4 shows  average  wave  number  spectra  for  each  of  two 
regions,  comprising  64  1/8  mesh  boxes,  i.e.,  a whole  mesh  box. 
Consequently,  the  FFT's  are  computed  on  8x8  arrays  and  wave 
number  spectra  are  produced  for  wave  numbers  up  to  six  (6<,4/2~ 
+.5<7).  Normalization  of  all  data  has  been  set  so  that  the 
average  value  is  the  d.c.  (zero  frequency)  component.  If 
classifications  were  to  proceed  from  this  point  using  these  data 
sets  for  training  purposes,  then  the  d.c.  component  would 
strongly  dominate  the  classification  procedure.  However,  if  we 
use  the  sample  standard  deviations  for  the  wave  number  spectra, 
computed  for  each  data  set  at  each  wave  number,  a different 
situation  emerges.  In  Figure  5,  we  have  plotted  the  ratios  of 
wave  number  spectra  to  standard  deviation.  Note  that  the 
relative  importance  of  the  d.c.  component  diminishes  greatly, 
particularly  for  the  case  representing  stratus  or  fog. 

Returning  to  Figure  4,  we  see  little  detail  in  either  wave 
number  spectrum  apart  from  the  disparity  between  the  d.c. 
component  and  components  at  wave  numbers  one  through  six.  Since 
VHR  data  were  not  available  to  allow  an  investigation  of  higher 
spatial  frequencies,  we  computed  two  FFT's,  one  for  the  64x64 
array  representing  the  stratus-stratocumulus  region,  and  the 
second  representing  the  stratus-fog  region;  these  data  are  shown 
in  Figure  6.  Note  that  wave  number  eight  in  this  figure 
corresponds  roughly  to  wave  number  one  in  Figures  4 and  5.  Exact 
correspondence  is  not  expected  because  additional  spectral 
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Figure  5.  Wave  number  spectra  from  two  whole  mesh  boxes 
normalized  by  class  conditional  standard 
deviations.  Each  spectrum  is  the  average  of 
64  8 x 8 spectra  from  HR  data. 
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Figure  6.  Wave  number  spectra  from  two  whole  mesh  boxes. 
Each  spectrum  Is  from  one  64  x 64  array. 
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components  from  different  orientations  and  frequencies  contribute 

to  the  average  for  the  64x64  based  wave  number  spectra.  To 

compute  the  spectra  of  Figure  6 and  scale  them  to  provide 

correspondence  with  the  spectra  of  Figure  4,  a multiplier  of 

eight  was  introduced.  This  brings  up  the  normalization  dichotomy 

mentioned  above.  Setting  the  d.c.  value  to  the  array  average  and 

using  the  same  factor  to  normalize  the  wave  number  spectra 

2 

ignores  the  factor  T which  should  also  be  applied  as  a 
multiplier  to  the  power  spectral  components.  Without  this  factor 
and  except  for  the  d.c.  which  remains  the  average  over  the 
window,  wave  number  spectra  will  decrease  as  the  number  of  terms 
in  the  FFT  increases.  However,  if  we  multiply  all  components  of 
the  power  spectrum  by  T , then  the  d.c.  component  will  increase 
as  the  number  of  terms  increases,  even  though  other  spectral 
components  may  fluctuate  about  the  same  general  level.  Rather 
than  consider  alternative  normalizations,  or  combine  the  two 
mentioned  above,  it  is  easily  seen  that  use  of  the  variance  data 
as  in  Figure  5 resolves  the  problem,  since  whatever  scale  factor 
is  used  will  appear  in  both  the  spectral  component  and  its 
standard  deviation  and  so  be  cancelled  in  the  ratio. 

While  there  are  obviously  more  details  in  the  spectra  of 
Figure  6 than  in  the  corresponding  spectra  of  Figure  4,  the 
former  are  not  exceedingly  more  useful  than  the  latter.  From  the 
low  wave  number  behavior  in  Figure  6,  it  appears  that  the 
stratus-stratocumulus  region  may  represent  a predominantly 
uniform  coverage,  whereas  the  flattening  of  the  stratus-fog 
spectrum  towards  low  wave  numbers  indicates  that  the  64x64  array 
(whole  mesh  box)  encompasses  several  correlation  intervals. 
Unfortunately,  absence  of  data  with  sufficient  resolution  ( VHR 
data)  makes  high  frequency  distinctions  mere  conjecture.  In  any 
spectral  analysis,  there  is  no  substitute  for  adequate  spectral 
and  spatial  resolution. 


-25- 


3 . 4 Summary  of  Results 


Since  the  completion  and  acceptance  of  our  R&D  Design  and 
Evaluation  Report  (Pickett  and  Blackman,  1976),  we  have  supported 
AFGWC  in  developing  the  capability  to  compute  wave  number  spectra 
for  cloud  classification.  Towards  this  end,  we  have:  (1)  aided 
the  FFT  implementation,  (2)  developed  a modified  wave  number 
spectrum  formulation,  (3)  recast  and  evaluated  certain  of 
Sikula's  data  (Sikula,  1973)  in  terms  of  this  modified  form,  (4) 
directed  the  processing  of  HR  data  from  display  tapes  to  test 
programs  and  demonstrate  feasibility,  and  (5)  established  an 
approach  to  be  used  in  the  development  of  a classifier  based  on 
wave  number  spectra.  In  this  last  regard,  computation  of  the 
FFT,  and  conversion  to  wave  number  spectra  produces  spectral 
features;  quadratic  discriminants  using  a diagonal  covariance 
matrix  will  classify  the  data  into  cloud  types.  Every  effort 
must  be  made  to  ensure  the  generation  of  a properly  labeled  truth 
set  of  images,  because  the  success  of  the  classifier  will  depend 
heavily  on  the  accuracy  with  which  class  conditional  means  and 
variances  are  estimated. 

4.  CONCLUSIONS  AND  RECOMMENDATIONS 

The  following  conclusions  and  recommendations  were  developed 
in  the  course  of  this  contract  and  are  discussed  in  detail  in  one 
or  the  other  of  the  two  reports: 

(1)  Further  studies  are  needed  to  establish  the  full  range  of 
need  for  imagery  processing  at  AFGWC,  particularly  with 
respect  to  automating  imagery  processing  now  done  visually. 
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